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A  computer  simulation  of  the  time-dependent 
Ginzburg-Landau  model  for  spinodai  decomposition 

Rolfe  Petschek  and  Horia  Metiu*’ 

Dtpanmtnt  of  Cktmistiy,  Unioartity  of  California,  Santa  Barbara.  California  93106 
{ Received  27  December  19(2;  accepted  I  March  1983) 

la  this  paper  we  discuss  the  results  of  computer  simulations  of  the  time-dependent  Ginzburg-Landau 
equation  appropriate  to  the  spinodai  decomposition  of  a  two-component  mixture  with  a  conserved  order 
parameter.  The  results  for  quantities  of  theoretical  interest,  such  as  the  probability  that  the  concentration  has 
a  particular  value  in  a  particular  ceil,  as  well  as  experimentally  accessible  quantities  such  as  the  structure 
taction,  are  presented.  For  reasons  of  computer  space,  the  simulation  was  done  in  two  dimensions.  The 
effects  of  incompressible  flow,  which  might  be  important  in  some  experimental  situations,  were  not  included 
in  these  calculations. 


f.  INTRODUCTION” 

The  term  spinodai  decomposition  is  used  to  describe 
the  process  by  which  a  two-component  mixture  in  a 
thermodynamically  unstable  region  separates  into  two 
equilibrium  phases.  •-*  This  is  a  process  of  technologi¬ 
cal  interest  since  it  takes  place  in  alloys,  *  polymers,* 
glasses, 1  and  binary  liquids.*  It  is  also  of  theoretical 
interest*-1*  since  it  provides  an  example  of  a  far  from 
equilibrium,  strongly  nonlinear  MdiffuSional“  process  in 
which  the  fluctuations  of  the  order  parameter  play  an 
important  role. 

Finally,  it  is  a  dynamic  process  taking  place  in  an 
unstable  phase,  below  the  critical  point  (we  have  in 
mind  a  phase  diagram  with  an  upper  critical  point), 
which  offers  an  opportunity  to  study  whether  the  dynam¬ 
ic  scaling  laws13  and  the  universality  behavior  observed 
above  T,  in  stable  phases  still  hold  below  T,  under  con¬ 
ditions  of  thermodynamic  instability. 

A  phenomenological  theory  for  the  most  probable 
evolution  of  the  concentration  during  decomposition  has 
been  proposed  by  <$h^i.*  This  was  Improved  by  Cook1* 
who  added  thermal  fluctuations .  Various  reformulations 
of  the  (^ffy-Cook  theory  were  presented  by  Langer,  u 
and  Metiu,  Kitahara,  and  Ross. u 

The  phenomenological  theory  is  heuristic  and  con¬ 
tains  a  number  of  approximations  whose  validity  has  not 
been  assessed  a  priori: 

(a)  It  is  assumed  that  the  probability  P[d,v),  t)  that 
the  concentration  at  a  point  r  has  the  value  e(r)  at  time 
t  satisfies  a  master  equation1,11'1*  with  the  transition 
rate 

wtc(r)-c<r)>6c<r|-aexp|-  J  drdr'Kij,  r')«c(r)  «c(r')j- 

.  x«xp{-^/Xc(r))-v8c(r))-F(c(r)))/2}.  (1) 

This  is  the  probability  that  the  concentration  changes  in 
unit  time  from  e(r)  to  c(r)  +  dc(r).  The  second  exponen¬ 
tial  in  Eq.  (1)  favors  changes  which  diminish  the  free 
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energy  F(c).  The  first  one  introduces  a  dynamic  con¬ 
straint  in  the  sense  that  the  probability  of  having  an  ex¬ 
cessively  large  change  8c(r)  per  unit  time  is  very  small, 
even  if  such  a  change  is  strongly  favored  by  the  free  en¬ 
ergy.  Even  though  a  theoretical  study  of  the  dynamic 
Ising  modelun>  showed  how  one  can  infer  a  transition 
rate  of  the  form  taken  in  Eq.  (1),  we  should  consider 
Eq.  (1)  to  be  untested  and  unproven  if  extended  to  other 
eases  such  as  spinodai  decomposition.  ' 

(b)  If  we  assume  that  r,  r')  is  a  number  A,  and  that 
A*1  is  a  small  parameter,  the  master  equation  can  be 
reduced1*11,1*  to  the  time  dependent  Ginzburg-Landau 
equation, 

»c(r,  t)/8f~MV*0F/S eir,  t )  ♦  tfr,  t) ,  (2) 

where  M  is  a  kinetic  coefficient  and  tj(r,  t )  is  a  Gaussian 
random  force  satisfying 

<ij(r,  t)  Hr,  /*)>  -  -  2*.  TM  6(t  - 1  ■)  v*  «(r  -  r*)  .  (3) 

The  latter  requirement  is  the  fluctuation  dissipation 
theorem.  Equations  (2)  and  (3)  are  the  Cahn-Cook  mod¬ 
el  for  spinodai  decomposition.  If  the  reduction  of  the 
master  equation  to  the  Langevin  equation  is  carried  out 
to  a  higher  order  in  A*1,  the  resulting  equation  is  rather 
different  from  Eq.  (3).  Furthermore,  if  A*1  is  not 
small  enough,  the  reduction  cannot  be  made  and  we  must 
work  with  the  master  equation.  There  is  so  far  no  ex¬ 
perimental  or  theoretical  proof  that  A’1  is  very  small; 
hence  that  Eq.  (2)  is  adequate. J)or-some-  models  of  in-  - 
terast.  i 

(e)  Finally,  all  theoretical  work  has  used  a  van  der 
Waals, 14  Chan-Hllllard, 11  Ginzburg-Landau,  “  Wilson11 
free  energy: 

F{r)  •  (*,  T/2)  f  dr  [/(e<r))  ♦  A(Vc(r))*]  -  F,  (4a) 

•<*»r/ 2)  f  dv {a c*(r) ♦  («/2) e4(r) «•  A(v c(r))*}  .  (4b) 

Here  kt  7f[c)/2  is  the  free  energy  density  of  a  uniform 
system  constrained1*  to  have  the  concentration  c. 
Throughout  this  paper  we  will  use  the  symbol  e  to  de¬ 
note  the  difference  between  the  concentration,  denoted 
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by  C,  and  the  critical  concentration"  and  will  refer  to 
both  loosely  as  "the  concentration. "  More  explicitly  e 
is  defined  as  that  variable  linear  in  the  actual  concen¬ 
tration  for  which  no  cubic  term  appears  in  the  function 
/.  The  gradient  term  takes  into  account  the  fact  that 
the  free  energy  is  affected  by  the  existence  of  a  spatial 
Inhomogeneity  in  the  concentration  field  and  that  the  free 
energy  at  a  given  point  depends  on  the  concentration  at 
neighboring  points.  The  gradient  term  stabilizes  inter¬ 
faces  between  phases  and  is  responsible  for  surface  ten¬ 
sion.  is  a  reference  free  energy  for  a  homogeneous 
system.  Equation  (4b)  is  obtained  from  Eq.  (4a)  by  ex¬ 
panding  /(c(r))  in  power  series  around  the  concentration 
of  the  homogeneous  reference  system. 


culations  to  the  experiment. 


Recent  work*3*11,1'  has  shown  that  Eq.  (4a),  and  the 
equation  obtained  from  (4a)  and  Eq.  (2),  are  a  satisfac¬ 
tory  starting  point  for  the  theory  of  static  and  some  dy¬ 
namic  critical  phenomena.  These  studies  considered, 
however,  only  steady-state  dynamic  phenomena  taking 
place  in  a  stable  region  of  the  phase  diagram,  above 
the  coexistence  line  and  close  to  the  critical  point. 
Therefore  their  success  does  not  imply  that  Eqs.  (2)- 
(4)  are  a  satisfactory  starting  point  for  a  description  of 
spinodal  decomposition  which  takes  place  in  the  unstable 
region  below  the  coexistence  line,  and  is  a  relaxation 
experiment  following  a  sudden  change  in  the  thermody¬ 
namic  state  which  can  be  carried  out  far  from  the  criti¬ 
cal  point. 

We  consider,  therefore,  that  the  assumptions  de¬ 
scribed  above  are  reasonable  but,  as  yet,  unproven. 
Attempts  to  Justify  them  are  hindered  by  the  fact  that 
the  competing  role  of  fluctuations,  nonlinear  interac¬ 
tions,  and  spatial  variation  of  the  order  parameter, 
make  these  equations  extremely  difficult  to  solve.  The 
most  ambitious  calculations  to  date  are  those  carried 
out  by  Langer,  Bar-on  and  Miller1  (LBM).  These  were 
extended  to  include  hydrodynamic  (mode-coupling)  ef¬ 
fects  by  Kawasaki  and  Ohta. 11  An  improved  LBM  type 
ansatz  was  recently  studied  by  Hill,  Metiu,  and 
Petschek.11  The  method  of  solution  involves  heuristic 
assumptions  whose  validity  is  very  difficult  to  assess. 

We  expect1  that  the  LBM  results  are  better  at  the  early 
times  of  the  decomposition  process.  Further  uncer¬ 
tainties  in  testing  these  equations  are  caused  by  the  pos¬ 
sibility11  that  the  existing  data  is  not  free  of  multiple 
scattering  and,  therefore,  It  should  not  be  compared  to 
single  scattering  calculations. 

Even  if  the  assumptions  above  are  correct  the  param¬ 
eters  which  go  into  the  theory  are  not  generally  easy  to 
measure  experimentally.  Thus  it  is  difficult  to  make 
critical  comparisons  between  theory  and  experiment. 

We  are  thus  in  a  situation  in  which  both  the  starting 
equations  and  their  solution  might  be  erroneous,  the  data 
are  not  yet  above  suspicion  and  comparisons  of  theory 
and  experiment  are  difficui..  For  these  reasons  we  have 
decided  that  some  progress  could  be  made  by  obtaining 
the  “exact"  solution  of  the  time -dependent  Ginzburg- 
Landau  equation  [Eqs.  (2),  (3),  and  (4b)]  by  a  novel  type 
of  computer  simulation.  We  can,  in  principle,  test  the 
existing  approximate  solutions  by  comparing  them  to 
our  “exact”  results,  and  we  can  also  compare  our  cal¬ 


The  simulation  consists  of  using  finite  differences, 
in  time  and  space,  to  turn  Eq.  (2)  into  an  atgebraie 
equation.  We  program  the  computer  to  generate  at  each 
“time  slice”  t,  and  spatial  point  r«  a  random  number 
»l(r„  t,)  taken  from  a  Gaussian  distribution  whose  width 
is  given  by  the  fluctuation-dissipation  theorem.  We 
then  solve  for  c(r„  tj  +  bt)  in  terms  of  these  random 
variables  ande(r„  tj).  We  then  repeat  this  procedure 
until  we  reach  some  latest  time.  A  complete  run  gives 
us  the  values  of  c(r„  tt)  at  all  space  points  r,  and  times 

tj.  We  call  the  set  {c(r„  /,),  i  *  1, 2 . a  *  0, 1, 2, ... } 

vith  a  “trajectory".  If  we  repeat  the  calculation  we  get 
second  trajectory,  which  differs  from  the  first  one  due 


to  stochastic  nature  of  the  quantity  q.  Running  many 
statistically  independent  trajectories  we  compute  the 
average  value  of  any  function  of  concentration. 

For  economic  reasons  we  have  carried  out  such  cal¬ 
culations  in  two  dimensions.  We  report  here  a  variety 
of  quantities:  the  evolution  of  the  concentration  in  time 
during  one  trajectory;  the  probability  P,(c,  t)  that  the 
concentration  at  a  given  point,  at  time  /„  equals  <r,  the 
probability  P,(e,,  c,\l)  that  at  time  t,  the  concentration 
at  site  i  is  c,  and  that  at  site  j  is  ct\  the  structure  factor 
£(&,/)  and  its  Fourier  transform  S(r,  /);  and  some  higher 
moments  of  P,(c,;/). 

Since  we  have  carried  out  two-dimensional  calcula¬ 
tions  we  cannot  compare  directly  with  three-dimensional 
experiments.  However  we  show  that  the  simulation  re¬ 
produces  the  qualitative  features  Observed  in  light  scat¬ 
tering1’4:  the  scattered  light  forms  a  ring  which  evolves 
in  time  so  that  its  intensity  grows  and  its  stze  contracts. 

The  present  theory1  predicts  that  if  the  starting  equa¬ 
tions  have  even  an  approximate  validity  for  physisorbed 
systems  on  smooth  surfaces  [chemisorbed  systems  can 
also  be  simulated  if  a  lattice-molecule  potential  is  added 
to  Eq.  (4)],  then  a  binary  adsorbed  system  will  undergo 
at  certain  concentrations  and  temperatures,  a  separa¬ 
tion  phase  transition  in  which  islands  of  high  concentra¬ 
tion  in  one  component  are  formed  in  a  sea  of  low  con¬ 
centration  of  the  same  component.  The  present  calcu¬ 
lation  describes  the  rate  of  island  formation.  While  no 
kinetic  measurements  have  been  made,  recent  work11 
has  shown  that  the  phase  diagram  of  a  Ar-Xe  mixture 
on  graphite  is  such  that  an  unmixing  process  such  as 
spinodal  decomposition  could  occur.  This  might  also 
happen  in  intercalated  compounds.11 

Finally  we  emphasize  that  the  present  computer  simu¬ 
lation  is  quite  different  from  prevtous  simulation  work 
which  used  the  Monte  Carlo  method11  to  study  the  de¬ 
composition  in  a  two-component,  nearest  neighbor,  ki¬ 
netic  Ising  model,  or  used  molecular  dynamics11  to 
study  the  decomposition  of  a  Lennard-Jones  fluid.  The 
connection  between  these  microscopic  theories  and  the 
Ginzburg- Landau  theory  used  by  the  phenomenological 
approach  is  unclear;  the  latter  must  be  a  coarse  grained 
(in  time  and  space)  version  of  the  former  models.  Fur¬ 
thermore,  in  the  specific  case  of  the  Ising  model,  it  is 
not  clear  whether  the  model  can  be  considered  a  repre¬ 
sentation  of  a  binary  fluid,  especially' far  away  from  the 


critical  temperature  when  universality  can  no  looser  be 
invoked.  In  the  case  o(  molecular  dynamics  the  time 
scale  is  of  the  order  of  picoseconds,  while  we  are  In¬ 
terested  in  hydrodynamic  time  scales.  Clearly  the 
simulation  presented  here  addresses  a  different  ques¬ 
tion  than  previous  work,  and  it  is  aimed  at  establishing 
with  accuracy  the  behavior  of  the  time  dependent  Ginz- 
burg- Landau  model. 

II.  THE  MOOEL  AND  THE  METHOD  OF  SOLUTION 
A.  The  model 

Equations  (2)  and  (4a)  lead  to 
»c(r,  atf -<*,  T/Z)ml  {»//»c(r,  t)  -  Kv*c{  r,  /)}  ♦  T,(r,  t) 

- TA/v*  { [“  flj+rtr.t)  .  (S) 

Here  c(r,  /)  is  the  difference  between  the  concentration 
C(r,  t)  at  the  point  r  at  time  t,  and  the  critical  concen¬ 
tration  Ct.  Another  important  parameter  la  the  con¬ 
centration  c0  of  the  system  before  the  quench.  If  c(*0 
(i.e.,  off  critical  quench),  it  is  better  to- work  with 
c'(r,  0  *  c(r,  t)  -  c#.  When  we  switch  to  this  new  vari¬ 
able,  a  quadratic  term  in  appears  in  Eq.  (5). 

Such  a  term  Is  thermodynamically  necessary  since,  for 
example,  if  c,<  0  the  quadratic  term  will  bias  the  sys¬ 
tem  towards  generating  a  larger  amount  of  the  equilib¬ 
rium  phase  having  lower  concentration.  If  c0»0,  the 
two  phases  are  generated  in  equal  amounts. 

If  K,  «,  and  i\  are  zero,  Eq.  (5)  reduces  to  the  ordi¬ 
nary  diffusion  equation.  Since  a*lhlf/ 8*0),.^.  a>  0  in 
the  thermodynamically  stable  region  of  the  phase  dia¬ 
gram  (i.e.,  above  the  coexistence  line)  and  the  "bars’* 
diffusion  coefficient  Da  *  Mbt  Ta  is  positive.  If  we  de¬ 
fine  the  spinodal  line  by  (»VTc)/8c*)-0,  then  the  bare 
diffusion  coefficient  becomes  negative  below  this  line. 
Equation  (S)  Implies  that  below  the  spinodal,  ll  K,u,r) 
are  zero,  a  concentration  fluctuation  grows  spontaneous¬ 
ly.  This  is  the  kinetic  manifestation  of  the  thermody¬ 
namic  instability. 

The  gradient  term  JCITcl*  is  introduced  into  the  the¬ 
ory  because,  as  the  system  becomes  spatially  inhomo¬ 
geneous,  some  free  energy  must  be  used  to  create 
boundaries  between  patches  having  different  concentra¬ 
tions.  The  gradient  term  accounts  for  the  appearance 
of  a  surface  tension  between  such  patches.  .This  term 
prevents  the  spontaneous  growth  of  the  high  vector 
Fourier  components  e »(f),  Since  such  components  re¬ 
sult  in  the  formation  of  sharp  domain  boundaries  with 
high  surface  tension. 

If  K*Q,  but  m,  n-0  (this  is  Cahn's  model*),  the  Fou¬ 
rier  transform  of  the  aquation  gives  it  *-MktT 
xkfta  +  IQ?)  e,(t)  and  the  Fourier  transform  of  the  effec¬ 
tive  diffusion  coefficient  is  D{k)  ■  Mb, T  k*(a  +  Kb*).  In 
the  unstable  region,  where  a  <  0,  the  Fourier  compo¬ 
nents  e^t)  with  *<V  \a\/K  grow  exponentially  while  the 
othe$deeay  exponentially. 

The  nonlinear  term  (it* 0)  becomes  mathematically 
important  as  the  decomposition  proceeds  and  c(r,  l ) 
grows.  Physically,  this  is  the  term  which  slows  down 


the  initial  exponential  growth  of  the  low  k  modes.  It 
also  contains  the  thermodynamic  information  concerning 
the  position,  on  the  phase  diagram,  of  the  two  final 
equilibrium  phases  towards  which  the  system  must 
evolve. 

The  random  force  ?)  simulates  thermal  fluctuations 
in  variables  not  Included  in  our  description.  Assuming 
that  they  are  Gaussian  is  justified  if  the  system  is  in 
local  thermodynamic  equilibrium  and  these  fluctuations 
relax  quickly. 

•  The  form  of  Eq.  (S)  is  determined  by  three  kinds  of 
assumptions. 

(a)  If  we  assume  that  A*1  is  small,  then  we  can  use 
Eq.  (2).  To  second  order  in  A*1  the  equation  for  »e/tt 
contains**  additional  terms  of  the  form  (3F/8c)*  and 
#*/'/  8c*,  which  are  ignored  here. 

(b)  The  free  energy  is  a  functional  of  c(r)  and  the 
gradient  terms  Te  •  Te,  v*e,  v4c,  etc. ,  which  are  ro- 
tationally  and  translationally  invariant.  Equation  (4)  is 
an  expansion  in  powers  of  c  and  Tc  to  the  lowest  physi¬ 
cally  reasonable  order.  All  the  terms  appearing  in  Eq. 
(4)  are  indispensable  for  a  qualitatively  correct  picture 
of  the  decomposition.  More  terms  might  be  necessary 
for  a  quantitatively  accurate  description. 

Above  the  critical  point  spatial  inhomogeneities  are 
created  by  thermal  fluctuations.  Since  the  phase  is 
thermodynamically  stable,  the  inhomogeneities  decay 
and  corresponding  concentration  fluctuations  on  large 
length  scales  are  small.  Furthermore,  near  the  critical 
point  the  "interfaces”  between  fluctuations  are  smooth 
sad  ITc  I*  tends  to  be  small.  This  implies  that  under 
these  conditions  Eq.  (4)  is  likely  to  be  quantiatively  ade¬ 
quate.  This  is  not,  however,  certain  in  spite  of  the  tact 
tiiat  the  model  has  been  so  successful  in  explaining  both 
static1*  and  dynamic11*1*  critical  phenomena.  These 
studies  have  concentrated  on  determining  the  critical 
exponents,  and  universality  assures  that  these  are  not 
sensitive  to  the  details  of  the  model. 

Spinodal  decomposition  takes  place  below  Tc  and  the 
spatial,  inhomogeneity  is  spontaneously  created  because 
the  system  Is  thermodynamically  unstable.  Thus  c(r,  t) 
grows  rapidly  and  higher  powers  of  c  might  be  necessary 
in  Eq.  (4).  If  the  quench  is  far  below  T,  there  is  no  a 
priori  reason  to  expect  that  (Tc)*  is  small;  higher  pow¬ 
ers  of  (Tc)*  might  be  needed  in  Eq.  (4)  as  well. 

Note  that  at  early  times  Eq.  (4),  and  the  resulting 
kinetic  Eq.  (5),  should  work  well  since  both  c  and  T*c 
are  expected  to  be  small.  For  this  reason  it  is  very 
important  to  carry  out  experiments  at  the  earliest  pos¬ 
sible  times  to  compare  to  the  theory  in  a  regime  in 
which  we  expect  it  to  work.  If  there  is  some  fundamen¬ 
tal  problem  with  these  equations,  more  basic  than  keep¬ 
ing  only  few  terms  in  the  expansion  mentioned  above,  it 
could  show  up  by  producing  disagreement  with  the  theory 
in  the  short.time  regime. 

(c)  A  third  type  of  approximation  is  made  by  assum¬ 
ing  that  the  concentration  is  the  only  relevant  field  in 
the  problem,  to  principle,  changes  in  concentration  are 
coupled  to  changes  in  energy  density  (or  temperature) 


and  Auid  velocity.  In  solids,  slsstle  strains  which  de¬ 
velop  during  decomposition  often  play  a  role. 

The  hydrodynamic  effects  hare  been  introduced  In  the 
LBMtheory  by  Kawsakl  and  Ohta.u  We  are  currently 
studying  such  effects  with  the  present  simulation  method 
and  hope  to  report  the  results  in  a  future  article.11  Hy¬ 
drodynamic  effects  are  not  expected  to  play  a  role  in 
alloys,  but  might  be  important  in  binary  liquids. 

The  l.mperature  field  has  two  possible  effects.  First 
the  temperature  and  concentration  fluctuations  are  cou¬ 
pled.  It  is  therefore  possible  that  the  temperature  fluc¬ 
tuations  will  affect  the  dynamic  behavior  of  the  system 
(provided  it  does  not  relax  too  quickly).  Second  if  the 
system  relaxes  adlabatically  the  mean  temperature  of 
the  system  will  change  during  the  decomposition.  This 
latter  effect  will  be  important  near  the  critical  point 
only  if  the  divergent  part  of  the  specific  heat  near  the 
critical  point  is  comparable  to  the  background  value  of 
the  specific  heat.  This  is  true  only  extremely  close  tp 
the  critical  point  in  typical  binary  fluids. u 


iteratively  for  c(r,  t  +  Af)  [the  previous  values,  p(r,  /), 
is  known].  In  practice,  between  two  and  four  Iterations 
were  usually  sufficient. 

The  generation  of  tj’  has  been  done  as  follows.  Using 
the  definition  Eq.  (9)  for  q'(r,  t )  we  compute 

<n'(r„  fi)j»'(n,  f*)>  -  -  2ft,  TM  A/  V*  «(r,  -  r,)  .  (10) 

Using  the  finite  difference  expressions  for  V*  [l.e., 

Eq.  (6)]  and  the  Kroneker  delta,  6ri,rj,  for  6(r,  -  rt), 
find 

•  <nWt)V(r,, /,)>•*  A*1  **™a/  (aal) 

if  r,»r,  and  /s 

<i?'(r„  /t)  ij'( r*.  /,)>  -  -  2  A'1  ft,  ra/Af  (lib) 

11  rt  and  rt  are  nearest  neighbors  and  *f*» 
and 

<n'(rl,fl)n'(rl,ft)>-0  (11c) 

otherwise. 


8.  The  method  of  solution 

The  time-dependent  Ginzburg- Landau  equation1  is  a 
nonlinear,  stochastic,  partial  differential  equation.  We 
have  approximated  it  numerically  in  the  following  fash¬ 
ion.  The  order  parameter  c(r)  is  given  at  a  discrete 
set  of  points  on  a  square  lattice.  This  corresponds  to 
taking  a  square  "Brillouin  zone”  rather  than  the  circu¬ 
lar  Brillouin  zone  favored  by  theoretical  treatments  of 
the  time-dependent  Ginzburg- Landau  equations.  This 
discrete  lattice  also  implies  that  the  underlying  point 
symmetry  is  only  a  square  symmetry  rather  than  the 
circular  symmetry  which  would  be  expected  in  many 
systems  of  physical  interest.  The  largest  such  lattice 
used  was  a  square  with  32  x32  sites,  due  to  limitations 
in  the  size  of  the  computer  available  for  this  work.  The 
Lapladan  was  replaced  by  its  finite  difference  exprvs- 

** - - - h'V^cc*^) 

( h***c(x,  yjw  -  4c(x,  y)  +c(x  ♦  A,  y)  *e(x  -  ft,  y) 


- +c(x,y+A)+c(x,  y-A) ,  {«) 

where  A  is  the  distance  between  the  centers  of  two 
neighboring  lattice  cells. 

This  replacement  reduces  Eq.  (5)  to  a  system  of  (32)* 
stochastic,  ordinary,  differential  equations  of  die  form 

£^.C({eir)})«.ij(r(,/),  (7) 


where  C  is  the  finite  difference  equivalent  of  the  right- 
hand  side  of  Eq.  (S)  and  {c(r)}«(c(r,), . . . ,  c(r„))  with 
m»  32  x32.  This  system  is  solved  by  writing 

[e(r,  /  ♦  At)  -  e(r,  f)y  Af 

.C({i(c(rf<+Af)*c(r,fj])*n'(r„f)  (8) 

with 


This  is  exactly  the  form  required  by  the  fluctuation- 
dissipation  theorem  applied  to  our  discrete  system. 

Practically  we  have  obtained  a  Gaussian  random  field 
with  the  covariance  given  by  Eqs.  (11)  by  defining 

n'(x,  y,  t)  -  A*1  l(t$x  *  ft,  y,  0  -  (£(*,  y, 

+  V*(x,  y + A,  /)  -  v*(x,  y,  <)] ;  (12) 

here  v,  and  v,  are  two  statistically  independent  Gaussian 
fields  with  the  covariance 

<vt(ru  /,)  v,(r„  /.)>  .  (i/,(rlP  /,)  v^tu  /t)>  -  2  A,  TM  (13a) 
if  rt  ■  r„  t|«t,  and 

<vi(r„  /,)  vilTf,  /*))  •  (t»j(rlt  /,)  vitr*,  /,)>  -  0  (13b) 

otherwise. 

Note  that  we  need  to  generate  one  value  of  ut  and  one 
for  vb  for  each  cell  of  the  square  lattice  [i.e.,  for  each 
discrete  value  of  the  pair  (x,  y)].  This  is  accomplished 
in  a  standard  way  by  using  a  packaged  generator  of  uni¬ 
formly  distributed  random  numbers.  In  our  calculations 
we  have  taken  the  initial  concentration  to  have  the  same 
value,  e*  in  every  cell.  This  corresponds  to  the  as¬ 
sumption  that  the  initial  temperature  is  so  high  that  there 
are  no  fluctuations  on  the  length  scales  of  interest.  IX 
the  initial  temperature  were  such  that  the  correlation 
length  is  comparable  to  our  lattice  spacing  ft  it  would 
be  necessary  to  have  fluctuations  in  c(r,  /«  0).  In  the 
appropriate  temperature  range  these  fluctuations  are 
approximately  Gaussian. 

We  will  assume  that  the  quench  is  instantaneous.  This 
can  be  experimentally  achieved  in  binary  liquids  by  using 
pressure  jumps.  We  are  currently  studying  the  behavior' 
of  the  system  under  the  influence  of  time  dependent  tem¬ 
perature  (or  pressure)  changes,  etc. 


*/»)«»  •  (» 

Using  the  prescription  described  below,  we  generate 
numerical  values  for  n^r,  /).  Then  we  solve  Eq.  (8) 


C.  The  quantities  computed 

There  are  a  number  of  quantities  that  give  some  in¬ 
sight  into  the  behavior  of  the  decomposing  system.  The 
ones  computed  here  are  discussed  below. 


The  "one  per  tide"  distribution  (unction  Prfc,  t)  is  the 
probability  that  the  concentration  in  one  cell,  at  time  l, 
has  the  value  c.  We  obtain  this  by  making  histograms 
of  the  values  o£  e  at  a  given  time  t,  for  a  sufficient  num¬ 
ber  of  trajectories.  In  particular  we  divide  the  range 
of  concentrations  into  which  the  concentration  falls  at 
any  given  time  into  a  number  of  evenly  spaced  intervals. 
Then  for  each  of  the  N  sites  on  the  lattice  and  the  K  tra¬ 
jectories  we  count  the  number  of  times,  M{e),  that  the 
concentration  falls  into  the  interval  centered  on  the  con¬ 
centration  c.  We  then  estimate  the  probability  as  Pt(c,  /) 

■  .l/(c)5t*‘  Ar*‘(Ac)"1,  where  Ac  is  the  width  of  the  con¬ 
centration  intervals. 

The  "two  particle’'  distribution  function  P,(c[,  is 
the  probability  that,  at  time  t,  the  concentration  in  cell 
i  is  c  j  and  the  concentration  in  cell  j  is  e\.  This  is  ob¬ 
tained  in  a  similar  way  by  making  histograms  for  the 
pair  (c|,  c‘j)  by  sampling  over  all  pairs  in  which  the  cells 
have  the  same  geometrical  position  with  respect  to  each 
other;  the  statistics  are  improved  by  including  in  the 
sampled  set  values  generated  by  many  trajectories. 

The  structure  factor  S(*,  0  is  defined  by 
S,(D  -  g  {Z  cf(0  c*(/)jo  /r‘  3i*‘  .  -  (15) 

Here 

c&)  -  £  expfik  •  rf) ,  (18) 

i 

where  r,  is  the  position  of  the  center  of  the  cell  i,  and 
k  is  the  two-dimensional  wave  vector,  constrained  to 
give  periodic  boundary  conditions.  That  is,  kt*(2v/L)nt 
and  ky  *(2v/L)ny,  where  n,,  ny  are  integers.  We  sum 
over  ail  the  vectors  k  having  the  length  k  and  over  all 
trajectories  and  divide  by  the  number  of  trajectories  91 
and  the  number  of  cells  S,  to  obtain  the  structure  fac¬ 
tor  per  cell. 

III.  RESUCTS~AND~DISCOSSION 
A.  The  choice  of  parameters 

Except  for  scale  factors,  the  equilibrium  behavior  of 
the  Ginzburg- Landau  equation  depends  on  two  dimen¬ 
sionless  parameters.  These  are  the  mean  field  corre¬ 
lation  length,  given  in  units  of  the  lattice  spacing  h,  £ur 

■  (  \a\/K)~l,*h?1  and  the  parameters  G - [(u/4) 
x(ff/2)"4/*]*/<4^,(2a)(*ar)a,  where  i  is  the  dimension¬ 
ality  of  the  space  and  a  *d- 1/2(4 -d).  The  paramet* 

G  is  the  Ginzburg  criterion  [see  Ref.  19(a),  Sec.  IU  6] 
which  is  used  to  estimate  how  important  the  fluctuations 
are  in  determining  the  values  of  the  equilibrium  quanti¬ 
ties  for  a  system  described  by  Ginzburg-Landau  free 
energy. 

For  the  simulations  discussed  in  this  paper  we  have 
chosen  values  {kt>2.5  (Ixr^O.S)  and  G>0.2.  Since  G 
is  small  we  expect  that  the  equilibrium  properties  of 
the  system  are  given,  to  good  approximation,  by  the 
mean  field  theory.  Moreover,  since  the  fluctuations 
are  not  dominant  we  do  not  expect  universal  behavior. 

To  establish  the  importance  of  the  nonlinear  effects 
during  the  spinodal  decomposition  we  compare  the  mean 


square  concentration  fluctuation  given  by  the  linear  the¬ 
ory  to 

C*,»[(Ct-C,)/2p,  (17) 

where  Ct  and  C,  are  the  concentrations  of  the  coexisting 
equilibrium  phases  at  the  temperature  at  which  the  de¬ 
composition  takes  place.  C„  gives  us  a  value  of  the 
concentration  fluctuation  at  which  we  know  that  nonlin¬ 
ear  effects  are  important.  The  mean  square  concentra¬ 
tion  fluctuations  for  the  linear  theory  are  given  by 

Cc(r,f)*)->/'rfkS(*,/) 

“ $;«(* » A;(VvO‘/*  exp(-  //4r0)  ,  (18) 

where  S^(A«0)  is  the  equilibrium  value  of  the  structure 
factor  at  the  temperature  T*  *  J*e+  I  ATI  (A7*  is  the  mag¬ 
nitude  of  the  temperature  quench  Tc-  T,  where  T  is  the 
temperature  at  which  the  decomposition  takes  place), 
A*«2'*'tv  ~4/lr(d/2)‘l  is  the  area  of  a  d-dimensional 
surface1,u’  divided  by  S*(/fe«  0),  f*  is  the  correlation 
length  at  the  temperature  7"  and 

T#»K/2Affl  (19) 

is  the  decomposition  time  scale  for  the  linear  theory 
(one-quarter  the  growth  rate  of  the  growing  mode).  The 
last  term  in  Eq.  (18)  is  the  asymptotic  (long  time)  value 
of  the  integral. 

Thus  we  see  that  the  length  of  time  during  which  the 
linear  theory  is  applicable  is  given  (in  units  of  r0)  by  the 
logarithm  of 

*-C^(rW/K,(*-0)/C,l.  (20) 

For  the  parameter  values  used  in  the  present  simula¬ 
tions,  the  value  of  R  is  10  (in  computing  this  we  used 
mean  field  values  for  £**  $•»,  C*»  and  C|.  As  the  ex¬ 
perimental3*  (and  theoretical)  value  of  R  in  three  di¬ 
mensions  is  16  we  expect  the  effect  of  the  nonlinearity 
in  the  present  simulation  to  be  comparable  to  that  in 
three  dimensions  near  the  critical  point. 

The  initial  concentration  in  all  the  simulations  pre¬ 
sented  here  was  taken  to  be  the  same  in  all  cells.  This 
Implies  that  the  system  is  initially  at  a  temperature  high 
above  the  critical  point  so  that  the  initial  fluctuations 
are  very  small  compared  to  those  occurring  during 
spinodal  decomposition.  In  one  simulation  the  initial 
concentration  was  taken  to  have  its  critical  value,  so 
that  e(t  ■  0)  ■  0.  In  the  other  we  use  c(f  * 0)  ■  -  0. 47 C„ 
(with  C„,  as  defined  in  Eq.  (18),  computed  within  mean 
field  Ginzburg-Landau  theory). 

B.  Results 

1.  The  concentration  pattern  and  its  growth 

Figures  1  and  2  show  qualitatively  the  concentration 
patterns  generated  by  one  trajectory,  at  three  different 
times  following  the  quench.  The  darkness  in  a  given 
cell  is  roughly  proportional  to  the  concentration  therein. 
The  darkest  cells  have  the  largest  positive  deviation  of 
the  concentration  from  the  critical  value,  and  the  light¬ 
est  ones  have  the  largest  negative  deviation.  In  the  off 
critical  quench  (Fig.  2)  the  cells  having  the  majority 


(c) 


phase  concentration  are  colored  lightly.  Note  that  as 
the  size  of  the  concentration  patches  increases  the 
boundary  between  them  becomes  sharper. 

2.  The  structure  factor  and  other  momenta 

In  Fig.  3  we  have  shown  typical  plots  of  the  dimen¬ 
sionless  structure  factor  per  unit  cell,  fC<le(*)l*)/W*f  F, 
where  tV  is  the  number  of  lattice  sites  and  K  is  defined 
in  Eq.  (4a).  The  structure  function  5  has  been  plotted 
as  a  function  of  the 

discrete  analog  of  the  magnitude  of  It.  In  particular, 

Ikl'  is  simply  the  square  root  of  the  Fourier  transform 
of  the  discrete  Laplacian  [Eq.  (6)].  It  has  been  made 
dimensionless  by  multiplying  by  the  lattice  spacing  h. 


FIG.  1.  This  figure  gives  a  qualitative  visual  impression  of 
the  concentration  distribution  during  spinodal  decomposition, 
in  critical  quench.  The  degree  of  darkness  is  proportional  to 
the  local  concentration.  The  light  (dark)  squares  are  placed 
where  the  local  concentration  is  smaller  (larger)  than  the  crit¬ 
ical  concentration.  The  Figs.  (a),  (b).  and  (c)  represent  the 
concentration  distribution  at  the  times  10  r0,  56  r,,  and  223  i#I 
respectively,  after  the  quench.  The  darkness  scale  has  been 
chosen  so  that  the  lightest  and  the  darkest  squares  represent 
places  where  c- 1. 1  C„,  that  have  the  maximum  concentration 
fluctuations.  The  mean  field  value  of  C„  (defined  in  Eq.  01)1 
was  used. 


It  la  essentially  identical  to  the  usual  definition  of  the 
magnitude  of  It  in  the  range  lkl'<2,  e.g.,  throughout 
the  range  of  interest.  The  structure  function  S  has  beer, 
averaged  over  all  wave  vectors  having  the  same  values 
of  Ikl'. 

The  error  bars  were  deduced  by  assuming  that  for 
each  trajectory  the  nonzero  or  imaginary  parts  of  the 
Fourier  transform  of  the  concentration  are  Gaussian 
random  variables.  This  is  the  expected  result  for  suf¬ 
ficiently  large  systems  since  each  Fourier  component 
of  the  order  parameter  is  the  sum  of  a  large  numbei  r>( 
essentially  independent  contributions  from  noninteract¬ 
ing  patches.  The  central  limit  theorem**  Implies  there¬ 
fore  that  the  sum  is  a  Gaussian  random  variable.  This 
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FIC.  2.  Same  as  Fig.  1  but  for  an  off-critical  quench.  The 
light  squares  are  used  for  the  majority  phase. 


argument,  of  course,  fails  for  integrals  over  the  Fouriet 
components  of  the  concentration,  e.g. ,  the  concentra¬ 
tion  at  a  given  point  is  not  necessarily  Gaussian.  We 
have  checked  the  above  assumption  by  calculating,  from 
the  results  of  our  simulation,  the  mean  square  fluctua¬ 
tion  for  the  square  of  each  Fourier  component.  We  find 
no  statistically  significant  deviations  from  Gaussian  be¬ 
havior. 

In  Figs.  4  and  5  we  have  shown  the  values  of  the 
structure  function  for  various  times  during  the  decor 
position.  In  these  figures  we  have  simply  drawn  smc 
curves  through  plots  similar  to  those  in  Fig.  3,  const s 


tent  with  the  error  bars.  For  these  times  we  have  not 
been  able  to  find  a  convincing  scaling  form  for  the  struc¬ 
ture  function  as  a  function  of  some  scaled  wave  vector. 

It  is  also  of  interest  to  Fourier  transform  the  struc¬ 
ture  factor  and  find  the  correlation  function  S(lrl) 

■  (c(r')c(r  +  r')).  This  function  is  plotted  in  Fig.  6  for 
several  times.  It  is  interesting  to  note  that  S(  I  rl)  is  a 
smooth  function  of  the  absolute  value  of  the  separation, 
for  sep  Motions  less  than  16  lattice  spacings  (this  is  half 
way  '  jss  the  finite  lattice).  It  therefore  appears  that 
th.  tape  of  the  underlying  square  lattice  is  relatively 
unimportant  to  the  structure  function  and  we  conjucture 


FIG.  3.  The  Fourier  transform  of  ths  structure  factor,  $<*,<)  as  a  function  of  *,  for  several  values  of  the  time  after  the  quench, 
(a)  the  critical  quench,  for  t «3»T,;,(b)  the  critical  quench,  for  1 » 223  t,j  (c)  the  off-critical  quench  for  !»5«  r|S  (d)  the  off- 
critical  quench;' for  t «  223  t,.  The  (one  standard  deviation)  error  bars  were  estimated  as  discussed  in  the  text. 


that  this  is  also  relatively  unimportant  to  the  other 
quantities  which  we  calculate.  It  Is  also  interesting  to 
note  that  the  dependence  of  the  structure  function  on 
Irl.  The  shape  at  small  Irl  is  easy  to  understand; 
nearby  points  tend  to  have  similar  values  of  the  concen¬ 
tration.  Because  of  the  conservation  law  there  is  a  sum 
rule  (fdrS(r)*Q)  and  therefore  S(l  rl)  must  be  negative 
at  larger  r.  However  we  find  a  persistent  peak  at  yet 
larger  r  which  contributes  appreciably  to  this  sum  rule. 
This  is  due  to  the  relative  sharpness  of  the  peak  in  the 
structure  function  as  a  function  of  *.  This  suggests  that 
some  care  should  be  taken  In  simulating  splnodal  de¬ 
composition  to  ensure  that  the  system  is  large  enough 
as  spatial  correlations  persist  appreciably  further  than 
the  naive  estimate  where  kmu  is  the  location  of 

the  maximum  in  the  structure  function  as  a  function  of  k 
at  that  time.  Note  that  for  the  latest  time  tor  which  we 
have  data  the  second  peak  In  the  structure  function  is 
aero  only  at  19  lattice  constants,  essentially  halfway 
across  our  finite  lattice. 

We  have  tested  the  validity  of  the  linear  theory  at 


short  times  by  comparing  the  results  of  our  simulation 
to  the  linear  equations  for  the  structure  factor.  We  find 
but  after  0. 5  and  1  linear  growth  times  (i.  e. ,  f/r0  ■  2 
or  4)  the  amplitude  of  the  fastest  growing  linear  mode 
already  exceeds  by  a  factor  of  2  that  given  by  the  full 
nonlinear  simulation. 

It  is  Interesting  to  ask  whether  or  not  this  disagree¬ 
ment  between  the  results  of  the  computer  simulation  and 
the  linear  theory  can  be  explained  simply  by  the  random 
phase  approximation  (RPA)  theory  for  the  interactions 
between  the  quickly  relaxing  high  k  modes  and  the  more 
slowly  relaxing  k  modes.  In  RPA  the  rate  of  change  of 
the  structure  function  Is  given  by® 

^-2ilf|k|*[(|k|^na,A)Sl*)  +  ll,  (21) 

where 

aarA*a+3/<5(r*0)  >  (22) 

Rather  than  solving  these  equations  and  comparing  the 
results  to  those  of  the  computer  simulation  we  have 


FIG.  5.  Same  as  Fig.  4,  but  for  the  off-crltlcal  quench. 


FIG .  4.  A  plot  of  the  structure  factor  5(*.  t)  •  K  <  I  e(k,  0 1 * 
XiVJbgD**  as  a  function  of  k.  The  value  of  *  is  computed  as 
described  in  Sec.  II B.  The  curves  were  visually  smoothed 
out,  consistent  with  the  estimated  error  bars,  (a)  Represents 
S(*.!>  for  a  critical  quench  and  various  times:  /  ■  10  rt,  a  solid 
line;  <>23  ra,  a  dotted  line;  i»38  r,,  dashed  Une;  ("36  r,t  a 
dot-dashed  Une.  (b)  Represents  S(h,f)  for  a  critical  quench 
for  various  times:  tv  77  r„  a  solid  Une;  i «  100  r„  a  dotted 
line;  f  *  127  r„  a  dashed  line,  (c)  Represents  Sik.t)  for  a 
critical  quench,  at  various  times:  t  •  136  r»,  a  dotted  Une;  t 
•  100  r#,  a  dotted  line;  f  ■  223  r,  a  dashed  line. 


chosen  to  examine  the  computer  simulation  to  see  when 
the  assumptions  of  the  RP  A  theory  breakdown.  It  seems 
reasonable  to  assert  that  this  breakdown  will  be  less  de¬ 
pendent  on  the  choice  of  initial  state  (which  is  not  easily 
controllable  in  real  experiments)  than  the  detailed  mean 
field  theory  dynamics. 


FIG.  6.  The  structure  factor  S(  I  r  I ,  /)  as  a  function  of  I  r  I , 
at  various  times  after  the  quench,  (a)  The  critical  quench: 
t«10  t,,  a  solid  line;  t-38  r„  a  dotted  line;  1-100  r„  a 
dashed  line;  1-223  T„  a  dot-dashed  line,  (b)  Same  as  (a),  but 
for  the  off-crlUcal  quench. 


fa  particular,  underlying  the  RPA  assumptions  is  that 
the  fluctuations  in  c(r)  are  Gaussian  (like  in  the  linear 
model,  but  with  the  width  renormalized  by  fluctuations). 
In  Fig.  7  we  give  plots  of  (c{r)*)  -  3  <c(r)*)*.  For  the 
off-critical  quench  we  have  also  given  <cJ);  for  a  critical 
quench  our  free  energy  functional  gives  <c])  >0.  If  the 
fluctuations  are  Gaussian  both  <c‘)  -  3  (c*)*  and  (c1)  are 
zero  and  the  non-Gaussian  nature  of  the  fluctuations  can 
be  deduced  from  the  dimensionless  ratios  (<c*)-3  (c*>*/ 
<c*>*  and  The  graphs  show  that  the  fluctua¬ 

tions  become  non-Gaussian  quite  early  in  the  decompo¬ 
sition,  indicating  that  nonlinear  effects  reach  quickly 
beyond  the  behavior  implied  by  the  random  phase  ap¬ 
proximation. 

It  is  difficult  to  assess  a  priori  the  statistical  uncer¬ 
tainty  in  the  one-point  moments  discussed  above  or  in 
various  other  quantities  which  are  discussed  below. 

This  is  because  the  statistical  distribution  of  such  quan¬ 
tities  at  a  single  point  is  a  complex  dynamical  quantity 
and  because  all  such  quantities  were  averaged  over  all 
the  sites  in  the  lattice.  Thus,  these  quantities  are  the 
sum  of  a  large  number  of  statistically  dependent  vari¬ 
ables  with  distributions  and  correlations  which  are  not 
known  a  priori.  Some  Idea  of  the  accuracy  of  these  mo¬ 
ments  can  be  obtained  by  calculating  the  value  of  mo¬ 
ments  such  as  <e(r)*),  which  are  required  by  symmetry 
.to  be  zero  in  the  critical  quench.  Such  tests  suggest 
that  the  statistical  accuracy  of  these  moments  and  those 
discussed  below  are  of  order  3%. 


3.  The  distribution  functions 

Another  quantity  of  interest  is  the  one  particle  dis¬ 
tribution  function  P,(c,  / )  defined  in  Sec.  AC.  This  is 
the  probability  that  the  concentration  at  an  arbitrary 
cell  has  the  value  c  at  time  t.  Numerical  results  for 
Pt(c,  t)  as  a  function  of  c,  at  various  times,  are  displayed 
in  Figs.  8  and  9.  The  initial  distribution  (/*  0)  has  non¬ 
zero  value  for  c«c0  only,  since  this  our  initial  condition. 
At  later  times  the  distribution  broadens  and  then  starts 
to  develop  two  peaks  near  the  mean-field  equilibrium 
values.  For  c0-0  (Fig.  8)  the  distribution  is  symme¬ 
tric  with  respect  to  0.  As  the  time  evolves  the  number 
of  cells  having  the  initial  concentration  is  substantially 
diminished  and  the  peak  values  become  more  and  more 
probable.  In  other  words  more  cells  tend  to  have  the 
concentration  near  the  mean-field  equilibrium  values. 

In  the  very  long  time  limit,  when  the  system  reaches 
equilibrium  we  expect  that  P^c;/-*)  consists  of  a  sum 
of  two  nearly  Gaussian  distributions  centered  around 


l  i  me 

FIG.  T.  We  show  here  the  dependence  of  time  of  various  quan¬ 
tities  related  to  the  RPA  approximation  (see  the  text).  (a)  For 
a  critical  quench  we  have:  Full  line:  We  plot  %  given  by  Eq. 
(22);  dotted  line:  We  plot  (  I  c(R,t)  I *)  /CL,  as  a  function  of 
time; dot-dashed  line.  We  plot  (c(R,f)J)/C},,  as  a  function  of 
time;  dashed  line:  We  plot  (<c(/t,<)*>  - 3  (c<ltft),)*>/C|,,  as  a 
function  of  time,  (b)  Same  as  (a),  but  tor  the  off-critlcal 
quench.  The  mean  field  value  of  C„,  defined  in  Eq.  (17),  was 
used. 
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FIG.  8.  The  probability  P(e,t)  of  finding  tho  concentration  c, 
la  «  coll  at  tho  time  t  after  a  critical  quench,  (a)  f « 10  r,j  (b) 
f-  38  r,;  (c)  ■  1  ■  223  r,.  See  the  text  for  a  dlacuaaloa  of  the 
deflaition  of 
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FIG.  9.  Same  aa  Fig.  8,  but  for  the  off-critical  quench. 


the  actual  equilibrium  concentrations,  each  Gaussian 
haring  a  width  given  by  !S,J(k)dk  where  S„(£)  is  the 
equilibrium  value  of  the  concentration  correlation  func¬ 
tion  in  the  corresponding  phase.  We  have  tried  to  fit 
Pt(c,  t)  at  earlier  times  with  such  a  double  Gaussian  form 
and  were  quantitatively  unsuccessful. 


FIG.  10.  Contour  plot*  o t  the  two- point  probability  distribution  P,(c,,ct{l)t  representing  the  probabiUty  that  at  Ume  t  the  con¬ 
centration  in  one  cell  la  ct  and  that  in  a  nearest  neighbor  Is  c,.  The  axes  represent  the  rallies  of  c,  and  e,.  The  tick  marks 
represent  the  mean  (laid  values  o£  C„.  The  solid  lines  are  lines  ol  constant  probabiUty  drawn  at  equally  spaced  Intervals 
CPle,,c,)-AP  a;  a  Integer).  Dotted  lines  have  been  drawn  at  Intermediate  intervals  [Pic,,  e,)«  API* +  0.5)1  In  regions  .where  the 
contours  are  widely  spaced.  The  value  of  &P  was  0.82  C^.  This  has  been  done  at  several  times;  (a)  ( *  10  t,j  (b)  l  -  38  r,;  (c) 
(-100  r,j  (d)  (-233  r,. 


Wa  have  also  computed  the  two-point  probability  dis¬ 
tribution  (unction  PAc„Cj-,t)  defined  in  Sec.  IIC.  This 
is  the  probabiUty  that,  at  time  t,  the  concentration  in 
the  cell  t  is  c,  and  that  in  cell  j  is  cr  In  Figs.  10  and 
11  we  show  contour  plots  joining  the  points  (c,,  c,)  of 
equal  probability.  This  should  be  read  as  a  topographic 
map  in  which  the  surface  of  the  "mountains"  is 
Pj(c„c,;/),  at  fixed/.  In  Figs.  10(a)-10(d)  we  present 
the  case  when  i  and  j  are  nearest  neighbors-. 

The  symmetry  of  these  figures  under  the  exchange 
c,~ct  is  a  consequence  of  the  definition  of  Pt.  The 


symmetry  under  c,—  -c„  ct - e,  exists  only  in  a 

critical  quench  for  a  free  energy  functional  which  is 
unchanged  by  e  —  -  c. 

Note  that  for  nearest  neighbors  the  probability  that 
one  cell  is  near  one  equilibrium  concentration  and  that 
a  neighboring  cell  is  at  the  other  equilibrium  concen¬ 
tration  is  relatively  small.  This  shows  that  the  incipient 
domain  walls  separating  the  areas  of  different  phases 
are  large  compared  to  lattice  spacing.  If  this  were  not 
the  case  we  would  expect  the  dynamics  to  depend  on  the 
underlying  lattice. 


FIG.  11.  In  this  figure  w*  hare  given  contour  plots  of  the  two- 
point  probability  distribution  for  a  number  of  pairs  of  cells  with 
different  separations  at  the  time  233  r,  (the  latest  time  shown 
in  Fig.  10).  A  description  of  these  plots  is  given  in  the  caption 
to  Fig.  10.  as  is  a  plot  for  nearest  neighbor  pairs  [Fig.  10(d)) 
at  this  time.  The  pairs  of  cells  in  these  plots  were  separated 
along  the  x  or  y  axes  of  the  computer  simulation  by  distances  of 
2  [Fig.  11(a)].  4  [Fig.  11(b)],  and  8  [Fig.  11(c)!  lattice  spacings. 


Note  also  the  formation  between  the  second  and  third 
times  of  a  minimum  in  die  distribution  function  near  c, 
*Cj*0. 

In  Figs.  ll(a)-ll(c)  we  have  shown  the  two-point 
probability  for  a  number  of  different  separation,  for  a 
critical  quench,  at  the  latest  time  we  have  simulated. 

Contour  plots  for  an  off  critical  quench  have  similar 
features.  However,  they  are  difficult  to  read  due  to  the 
very  large  probability  that  both  points  have  concentra¬ 
tions  near  that  of  the  majority  phase  and  are  not  re-  _ 
produced  here.  (j 

It  is  important*0’"  to  have  a  good  description  of  this 
two-point  correlation  function  in  terms  of  a  small  num-  i*' 
ber  of  parameters.  A  comparison  of  various  parame¬ 
terization!  which  have  been  used  in  theoretical  work 
with  the  actual  measured  distributions  is  currently  in 
progress. 


ACKNOWLEDGMENTS 

This  work  has  been  supported  in  part  by  die  Office  of 
Naval  Research.  We  are  grateful  to  Doug  Scalapino, 
Susan  Hill,  Paul  Meakin,  Jim  Langer,  and  David  Cannell 
for  useful  discussions. 


*J.  Cahn,  Metall.  Soc.  A.I.M.E.  242,  166  (1968). 

L  S.  Longer,  in  Fluctuations,  Instabilities  and  Phase  Tran- 
^  silions,  edited  by  T.  Rlste  (Plenum,  New  York,  1975). 

*H.  Metiu,  K.  Kltahara,  and  J.  Ross,  In  Fluctuation  phenom¬ 
ena,  edited  by  E.  W.  Montroll  and  J.  L.  Lebowitx  (North- 
Hoiland,  New  York,  1979). 

4F.  F.  Abraham,  Phya.  Rep  53,  95  (1979). 

*K.  B.  Rundman  and  J.  E.  Hilliard,  Acta  Metall.  15,  1025 
(196?h  S.  Agarwal  and  H.  Herman,  Scr.  Metall.  7,  503 
(197 3)i  R.  AlcuBa  and  A.  Bonfiglio,  Acta  Metall.  22,  399 
(1974)}  V.  Gerald  and  \V.  Marts,  Scr.  Metall.  1,  33  (1967); 


R.  Petschek  and  H.  Metiu:  Ginzburg- Landau  moedl  for  spinodal  decomposition 


T.  L.  Bartel  and  K.  B.  Rundman,  Metals.  Tran*,  a,  1887 
(1975);  Y.  L\  Hilliard  and  K.  B.  Rundmar.  Scr.  1. 

37  (1967). 

*G.  T.  Peke  and  W.  Prias,  Macromoleculcs  7.  537  U9.'4>. 
tG.  P.  Kcilsen,  phya.  Chem.  Glasses.  10,  $4  (I960);  M. 
Tomozawa,  R.  K.  McCronc,  and  H.  Herman,  J.  Am.  Ceram. 
Soc.  55.  52  (1970);  N.  S.  Andreev,  G.  G.  Boiko,  and  K.  A. 
Bokov,  Noncryst.  Solids  5,  41  (1970). 

*J.  S.  Huang,  W.  I.  Goldburg,  and  A.  W.  Bjerkaas,  Phya. 

Rev.  Lett.  39,  931  (1974);  A.  J.  Schwartz,  J.  S.  Huang,  and 
U'.  I.  Goldburg,  J.  Chem.  Phys.  63,  1S47  (1975);  N.  Wong 
andC.  M.  Knobler,  MiW.  6*-.  4707  (1977). 

*J.  W.  Cahn.  Acta.  Mctall.  9.  795  (1961);  10,  179  (1962). 

**H.  E.  Cook.  Acta.  Metall.  IS.  297  (1970). 

“j.  S.  Longer,  Ann.  Phys.  M.Y.  54.  >58  (1969);  65.  S3 
(1971). 

**(a)  H.  Metiu,  K.  Kitahara,  and  J.  Ross,  J.  Chem.  Phys.  65, 
393  (1978);  (b)  63,  5116  (1975). 

UP.  C.  Hohenberg  and  B.  I.  Halperln,  Rev.  Mod.  Phys.  43, 

433  (1977). 

UJ.  D.  van  der  Waals,  Z.  Phys.  Chem.  (Leipzig)  25,  258 
(1894). 

UJ.  W.  Cahn  and  J.  E.  Hilliard,  J.  Chem.  Phys.  28.  258 
11938). 

‘*V.  L.  Ginzburg  and  L.  D.  Landau,  Zh.  Eksp.  Teor.  FI*.  20. 
1064  (1950). 

>tK.  G.  Wilson.  Phys.  Rev.  Lott.  28  ,  580  (1972);  K.  G. 

Wilson  and  J.  Kogut,  Phys.  Rep.  12,  75  (1974). 

UJ.  S.  Longer,  Physics  (Utrecht)  73,  61  (1974). 

>f(a)  S.  K.  Ma,  Modern  Theory  of  Critical  Phenomena  (Benja- 
mon,  Reading,  Mass,  1976).  (b)  D.  4.  Amit,  Field  Theory, 
The  Rcnormatisaiioa  Croup  and  Critical  Phenomena  (McGraw- 
Hill.  New  York,  1978). 

**4.  S.  Longer,  M.  Bar-on,  and  H.  O.  Miller,  Phys.  Rev.  A 


11,  1417  (1975). 

:iK.  Kawasaki  and  T.  Ohta,  Prog.  Teor.  Phys.  59,  3C2  (1978). 

»S.  C.  Hill,  H.  Metiu  and  R.  G.  Petschek  (unpublished). 

nlt  is  hoped  that  ongoing  experiments  carried  out  by  D.  S. 

,  Cannoil  and  C.  M.  Jefferson  with  a  mixture  3-mcthylpentane- 
nitroethane  in  which  the  refractive  indices  of  the  components 
are  nearly  equal,  will  provide  data  free  of  multiple  scattering. 

,4W.  I.  Goldburg  and  J.  S.  Huang,  in  Fluctuations ,  tnstabilites 
and  Phase  Transitions,  edited  by  T.  Riste  (Plenum,  New 
York.  1975). 

Bj.  Bohr,  M.  Nielson,  J.  P.  McTague.  J.  Als-Nlclsen,  and 
K.  Kjaer  (preprint). 

MPor  a  phase  diagram  see  J.  A.  Safran,  Phys.  Rev.  Lett.  44, 
937  (I960);  Dr.  Safran  has  suggested  that  spinodal  decomposi¬ 
tion  might  take  place  in  intercalated  compounds  (private 
communication). 

nA.  B.  Bortz,  51.  H.  Kalos,  J.  L.  Lcbowitz,  and  M.  A. 
Zendejas,  Phys.  Rev.  B  10,  535  (1974);  J.  Marro,  A.  B. 
Bortz,  M.  H.  Kalos,  andj.  L.  Lcbowitz,  ibid.  12,  2000 
(1975);  M.  Rao,  M.  11.  Kalos,  J.  Marro,  and  J.  L. 

Lcbowitz,  ibid.  13,  4328  (1976);  A.  Sur.  J.  L.  Lcbowitz.  J. 
Marro,  and  M.  H.  Kalos.  ibid.  IS,  3014  (1977). 

**M.  R.  Mruzik,  F.  P.  Abraham,  and  G.  M.  Pound,  J.  Chem. 
Phys.  69,  3462  (1978). 

**[t  is  Implicitly  assumed  in  Eq.  (4a)  that  the  coexistance  line 
is  symmetric  with  rcsiiect  to  T^.  This  is  generally  true  (or 
concentration  close  to  ct.  This  assumption  can  be  easily 
removed. 

**H.  Metiu  (unpublished). 

]|R.  Petschek  and  H.  Metiu  (to  be  published). 

ttH.  Guttinger  and  D.  S.  Connell.  Phys.  Rev.  A  6,  3188 
(1961). 

aA.  1.  Khlnchin,  Mathematical  Foundation  of  Statistical  Me¬ 
chanics  (Dover,  New  York,  1949),  Chap.  V. 


